Plots for the dataset of DMTA on Devionan sharks
Ivan Calandra
2023-12-12 15:35:45 CET
1 Goal of the script
The script plots all SSFA variables for the Devonian shark dataset.
Input Rbin data file must be located in
“./analysis/derived_data”.
Plots will be saved in “./analysis/plots”.
The knit directory for this script is the project directory.
2 Load packages
pack_to_load <- sort(c("R.utils", "tidyverse", "factoextra", "ggplot2", "ggpubr", "RColorBrewer",
"knitr", "rmarkdown", "grateful"))
sapply(pack_to_load, library, character.only = TRUE, logical.return = TRUE) factoextra ggplot2 ggpubr grateful knitr R.utils
TRUE TRUE TRUE TRUE TRUE TRUE
RColorBrewer rmarkdown tidyverse
TRUE TRUE TRUE 3 Read in data
3.1 Get name and path of input file
[1] "analysis/derived_data/DMTA-Ctenacanths.Rbin"3.2 Read in Rbin file
'data.frame': 160 obs. of 43 variables:
$ Specimen : chr "CC" "CC" "CC" "CC" ...
$ Tooth : chr "A" "A" "A" "A" ...
$ Location : chr "loc1" "loc1" "loc1" "loc1" ...
$ Objective : chr "100x" "100x" "100x" "20x" ...
$ Measurement : chr "meas1" "meas2" "meas3" "meas1" ...
$ Position : chr "top" "top" "top" "top" ...
$ NMP : num 3.03 3.05 3.34 12.54 11.63 ...
$ NMP_cat : Ord.factor w/ 3 levels "<10%"<"10-20%"<..: 1 1 1 2 2 2 1 1 1 3 ...
$ Sq : num 1.22 1.18 1.16 7.28 7.29 ...
$ Ssk : num -0.332 -0.345 -0.355 0.83 0.814 ...
$ Sku : num 2.53 2.45 2.42 3.15 3.14 ...
$ Sp : num 2.74 2.73 2.69 22.04 22.01 ...
$ Sv : num 3.42 3.25 3.17 17.88 19.03 ...
$ Sz : num 6.16 5.98 5.86 39.92 41.04 ...
$ Sa : num 0.995 0.969 0.956 5.8 5.803 ...
$ Smr : num 5.865 4.991 4.781 0.189 0.183 ...
$ Smc : num 1.5 1.47 1.43 11.42 11.45 ...
$ Sxp : num 2.62 2.5 2.48 8.72 8.84 ...
$ Sal : num 20.7 20.4 20.3 91.8 91.5 ...
$ Str : num 0.557 0.53 0.537 NA NA ...
$ Std : num 58.5 58.5 58.5 58.5 58.7 ...
$ Ssw : num 0.658 0.658 0.658 3.295 3.295 ...
$ Sdq : num 0.179 0.176 0.176 0.252 0.257 ...
$ Sdr : num 1.51 1.47 1.46 2.97 3.08 ...
$ Vm : num 0.0382 0.033 0.032 0.4189 0.414 ...
$ Vv : num 1.54 1.5 1.47 11.84 11.86 ...
$ Vmp : num 0.0382 0.033 0.032 0.4189 0.414 ...
$ Vmc : num 1.21 1.18 1.15 5.85 5.88 ...
$ Vvc : num 1.39 1.35 1.32 11.35 11.37 ...
$ Vvv : num 0.148 0.148 0.147 0.487 0.496 ...
$ Maximum.depth.of.furrows: num 2.33 2.29 2.29 41.48 39.05 ...
$ Mean.depth.of.furrows : num 0.834 0.808 0.792 18.822 18.342 ...
$ Mean.density.of.furrows : num 2027 2022 2034 1354 1319 ...
$ First.direction : num 37 37 37 90 90 ...
$ Second.direction : num 56.5 90 90 37 37 ...
$ Third.direction : num 90 56.6 56.6 56.5 56.5 ...
$ Texture.isotropy : num 47.3 51.9 58.7 21 21 ...
$ epLsar : num 0.00119 0.00113 0.00124 0.00347 0.00337 ...
$ NewEplsar : num 0.0172 0.0172 0.0173 0.0163 0.0164 ...
$ Asfc : num 1.85 1.79 1.82 6.64 6.88 ...
$ Smfc : num 86.2 48.8 52 12.9 12.9 ...
$ HAsfc9 : num 0.249 0.233 0.172 0.364 0.363 ...
$ HAsfc81 : num 0.591 0.582 0.516 0.486 0.459 ...
- attr(*, "comment")= Named chr [1:36] "%" "µm" "<no unit>" "<no unit>" ...
..- attr(*, "names")= chr [1:36] "NMP" "Sq" "Ssk" "Sku" ...4 Exclude surfaces with NMP ≥ 20%
Surfaces with more than 20% NMP should not be analyzed.
sharks_nmp20 <- filter(sharks, NMP_cat != "≥20%")
sharks_nmp20$NMP_cat <- factor(sharks_nmp20$NMP_cat)
str(sharks_nmp20)'data.frame': 144 obs. of 43 variables:
$ Specimen : chr "CC" "CC" "CC" "CC" ...
$ Tooth : chr "A" "A" "A" "A" ...
$ Location : chr "loc1" "loc1" "loc1" "loc1" ...
$ Objective : chr "100x" "100x" "100x" "20x" ...
$ Measurement : chr "meas1" "meas2" "meas3" "meas1" ...
$ Position : chr "top" "top" "top" "top" ...
$ NMP : num 3.03 3.05 3.34 12.54 11.63 ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 1 1 1 2 2 2 1 1 1 2 ...
$ Sq : num 1.22 1.18 1.16 7.28 7.29 ...
$ Ssk : num -0.332 -0.345 -0.355 0.83 0.814 ...
$ Sku : num 2.53 2.45 2.42 3.15 3.14 ...
$ Sp : num 2.74 2.73 2.69 22.04 22.01 ...
$ Sv : num 3.42 3.25 3.17 17.88 19.03 ...
$ Sz : num 6.16 5.98 5.86 39.92 41.04 ...
$ Sa : num 0.995 0.969 0.956 5.8 5.803 ...
$ Smr : num 5.865 4.991 4.781 0.189 0.183 ...
$ Smc : num 1.5 1.47 1.43 11.42 11.45 ...
$ Sxp : num 2.62 2.5 2.48 8.72 8.84 ...
$ Sal : num 20.7 20.4 20.3 91.8 91.5 ...
$ Str : num 0.557 0.53 0.537 NA NA ...
$ Std : num 58.5 58.5 58.5 58.5 58.7 ...
$ Ssw : num 0.658 0.658 0.658 3.295 3.295 ...
$ Sdq : num 0.179 0.176 0.176 0.252 0.257 ...
$ Sdr : num 1.51 1.47 1.46 2.97 3.08 ...
$ Vm : num 0.0382 0.033 0.032 0.4189 0.414 ...
$ Vv : num 1.54 1.5 1.47 11.84 11.86 ...
$ Vmp : num 0.0382 0.033 0.032 0.4189 0.414 ...
$ Vmc : num 1.21 1.18 1.15 5.85 5.88 ...
$ Vvc : num 1.39 1.35 1.32 11.35 11.37 ...
$ Vvv : num 0.148 0.148 0.147 0.487 0.496 ...
$ Maximum.depth.of.furrows: num 2.33 2.29 2.29 41.48 39.05 ...
$ Mean.depth.of.furrows : num 0.834 0.808 0.792 18.822 18.342 ...
$ Mean.density.of.furrows : num 2027 2022 2034 1354 1319 ...
$ First.direction : num 37 37 37 90 90 ...
$ Second.direction : num 56.5 90 90 37 37 ...
$ Third.direction : num 90 56.6 56.6 56.5 56.5 ...
$ Texture.isotropy : num 47.3 51.9 58.7 21 21 ...
$ epLsar : num 0.00119 0.00113 0.00124 0.00347 0.00337 ...
$ NewEplsar : num 0.0172 0.0172 0.0173 0.0163 0.0164 ...
$ Asfc : num 1.85 1.79 1.82 6.64 6.88 ...
$ Smfc : num 86.2 48.8 52 12.9 12.9 ...
$ HAsfc9 : num 0.249 0.233 0.172 0.364 0.363 ...
$ HAsfc81 : num 0.591 0.582 0.516 0.486 0.459 ...
- attr(*, "comment")= Named chr [1:36] "%" "µm" "<no unit>" "<no unit>" ...
..- attr(*, "names")= chr [1:36] "NMP" "Sq" "Ssk" "Sku" ...5 Separate data acquired with different objectives
#20x objective
sharks_nmp20_20x <- filter(sharks_nmp20, Objective == "20x")
sharks_nmp20_20x$Objective <- factor(sharks_nmp20_20x$Objective)
str(sharks_nmp20_20x)'data.frame': 69 obs. of 43 variables:
$ Specimen : chr "CC" "CC" "CC" "CC" ...
$ Tooth : chr "A" "A" "A" "A" ...
$ Location : chr "loc1" "loc1" "loc1" "loc3" ...
$ Objective : Factor w/ 1 level "20x": 1 1 1 1 1 1 1 1 1 1 ...
$ Measurement : chr "meas1" "meas2" "meas3" "meas1" ...
$ Position : chr "top" "top" "top" "bottom" ...
$ NMP : num 12.5 11.6 11.5 13.6 13.5 ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 2 2 2 2 2 2 1 1 1 1 ...
$ Sq : num 7.28 7.29 7.3 7.91 7.89 ...
$ Ssk : num 0.83 0.814 0.794 -0.714 -0.717 ...
$ Sku : num 3.15 3.14 3.13 4.11 4.13 ...
$ Sp : num 22 22 21.9 17.5 17.4 ...
$ Sv : num 17.9 19 19.6 33.3 33.4 ...
$ Sz : num 39.9 41 41.5 50.7 50.9 ...
$ Sa : num 5.8 5.8 5.81 6.09 6.07 ...
$ Smr : num 0.189 0.183 0.179 0.319 0.303 ...
$ Smc : num 11.42 11.45 11.42 9.81 9.78 ...
$ Sxp : num 8.72 8.84 8.94 19.26 19.17 ...
$ Sal : num 91.8 91.5 91.3 66.3 66.1 ...
$ Str : num NA NA NA NA NA ...
$ Std : num 58.51 58.75 58.51 3.74 3.75 ...
$ Ssw : num 3.29 3.29 3.29 4.39 4.39 ...
$ Sdq : num 0.252 0.257 0.257 0.33 0.33 ...
$ Sdr : num 2.97 3.08 3.09 4.98 4.97 ...
$ Vm : num 0.419 0.414 0.413 0.231 0.231 ...
$ Vv : num 11.8 11.9 11.8 10 10 ...
$ Vmp : num 0.419 0.414 0.413 0.231 0.231 ...
$ Vmc : num 5.85 5.88 5.89 6.49 6.46 ...
$ Vvc : num 11.35 11.37 11.33 8.82 8.8 ...
$ Vvv : num 0.487 0.496 0.507 1.218 1.216 ...
$ Maximum.depth.of.furrows: num 41.5 39.1 36.9 47.2 47 ...
$ Mean.depth.of.furrows : num 18.8 18.3 18.3 16.9 16.9 ...
$ Mean.density.of.furrows : num 1354 1319 1322 1443 1442 ...
$ First.direction : num 90 90 90 143 143 ...
$ Second.direction : num 37.00583 37.00101 37.00585 0.00521 0.0164 ...
$ Third.direction : num 56.5 56.5 56.5 20.8 20.8 ...
$ Texture.isotropy : num 20.96 20.96 20.92 9.82 9.8 ...
$ epLsar : num 0.00347 0.00337 0.0033 0.00213 0.00211 ...
$ NewEplsar : num 0.0163 0.0164 0.0164 0.0181 0.018 ...
$ Asfc : num 6.64 6.88 7.01 9.35 9.23 ...
$ Smfc : num 12.9 12.9 12.9 12 12 ...
$ HAsfc9 : num 0.364 0.363 0.37 0.265 0.276 ...
$ HAsfc81 : num 0.486 0.459 0.462 0.368 0.364 ...
- attr(*, "comment")= Named chr [1:36] "%" "µm" "<no unit>" "<no unit>" ...
..- attr(*, "names")= chr [1:36] "NMP" "Sq" "Ssk" "Sku" ...#100x objective
sharks_nmp20_100x <- filter(sharks_nmp20, Objective == "100x")
sharks_nmp20_100x$Objective <- factor(sharks_nmp20_100x$Objective)
str(sharks_nmp20_100x)'data.frame': 75 obs. of 43 variables:
$ Specimen : chr "CC" "CC" "CC" "CC" ...
$ Tooth : chr "A" "A" "A" "A" ...
$ Location : chr "loc1" "loc1" "loc1" "loc2" ...
$ Objective : Factor w/ 1 level "100x": 1 1 1 1 1 1 1 1 1 1 ...
$ Measurement : chr "meas1" "meas2" "meas3" "meas1" ...
$ Position : chr "top" "top" "top" "bottom" ...
$ NMP : num 3.03 3.05 3.34 9.84 9.78 ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 1 1 1 1 1 1 2 2 2 1 ...
$ Sq : num 1.22 1.18 1.16 1.66 1.66 ...
$ Ssk : num -0.332 -0.345 -0.355 1.53 1.454 ...
$ Sku : num 2.53 2.45 2.42 8.22 7.72 ...
$ Sp : num 2.74 2.73 2.69 8.77 8.39 ...
$ Sv : num 3.42 3.25 3.17 3.39 3.45 ...
$ Sz : num 6.16 5.98 5.86 12.16 11.84 ...
$ Sa : num 0.995 0.969 0.956 1.14 1.154 ...
$ Smr : num 5.865 4.991 4.781 0.361 0.405 ...
$ Smc : num 1.5 1.47 1.43 1.3 1.32 ...
$ Sxp : num 2.62 2.5 2.48 2.77 2.76 ...
$ Sal : num 20.7 20.4 20.3 15.3 15.2 ...
$ Str : num 0.557 0.53 0.537 0.408 0.396 ...
$ Std : num 58.5 58.5 58.5 84.5 84.5 ...
$ Ssw : num 0.658 0.658 0.658 0.658 0.658 ...
$ Sdq : num 0.179 0.176 0.176 0.409 0.395 ...
$ Sdr : num 1.51 1.47 1.46 4.86 4.71 ...
$ Vm : num 0.0382 0.033 0.032 0.2068 0.2047 ...
$ Vv : num 1.54 1.5 1.47 1.51 1.53 ...
$ Vmp : num 0.0382 0.033 0.032 0.2068 0.2047 ...
$ Vmc : num 1.21 1.18 1.15 1.16 1.18 ...
$ Vvc : num 1.39 1.35 1.32 1.35 1.37 ...
$ Vvv : num 0.148 0.148 0.147 0.154 0.155 ...
$ Maximum.depth.of.furrows: num 2.33 2.29 2.29 4.51 4.52 ...
$ Mean.depth.of.furrows : num 0.834 0.808 0.792 1.292 1.286 ...
$ Mean.density.of.furrows : num 2027 2022 2034 2211 2229 ...
$ First.direction : num 37 37 37 90 90 ...
$ Second.direction : num 56.5 90 90 84.3 84.3 ...
$ Third.direction : num 90 56.6 56.6 78.7 78.7 ...
$ Texture.isotropy : num 47.3 51.9 58.7 72.8 68 ...
$ epLsar : num 0.00119 0.00113 0.00124 0.00152 0.0018 ...
$ NewEplsar : num 0.0172 0.0172 0.0173 0.0176 0.0177 ...
$ Asfc : num 1.85 1.79 1.82 9.4 7.54 ...
$ Smfc : num 86.2 48.8 52 37.9 62.8 ...
$ HAsfc9 : num 0.249 0.233 0.172 3.11 1.59 ...
$ HAsfc81 : num 0.591 0.582 0.516 4.706 5.056 ...
- attr(*, "comment")= Named chr [1:36] "%" "µm" "<no unit>" "<no unit>" ...
..- attr(*, "names")= chr [1:36] "NMP" "Sq" "Ssk" "Sku" ...6 Plot each surface parameter in a boxplot
6.1 Define variables
Here we define which columns are used for the boxplots.
# Columns to be used to group on the x-axis
x_sp <- "Specimen"
x_pos <- "Position"
# Columns to be used on the y-axis
y_ISO_height <- colnames(sharks)[9:15]
y_ISO_vol <- colnames(sharks)[25:30]
y_ISO_others <- colnames(sharks)[16:24]
y_furrow_diriso <- colnames(sharks)[31:37]
y_SSFA <- colnames(sharks)[38:43]
# colors
grp_colors <- "Tooth"
# shapes
grp_shapes <- "NMP_cat"6.2 Plotting function
custom_boxplot <- function(dat, x_axis, y_axis = "Value",
group_col = grp_colors, group_shape = grp_shapes,
facet_var = "Parameter", plot_title){
# Define aesthetics
p_out <- ggplot(dat, aes(x = .data[[x_axis]], y = .data[[y_axis]])) +
# Boxplots:
# hide outliers (all points are shown with geom_point() below)
geom_boxplot(outlier.shape = NA) +
# Points:
# Add layers of shapes and colors for points
# Jitter points
geom_point(mapping = aes(shape = .data[[group_shape]], color = .data[[group_col]]),
position = "jitter", size = 3, alpha = 0.7) +
# Wrap around parameters (i.e. 1 subplot per parameter) with free y-scales
facet_wrap(as.formula(paste0("~", facet_var)), scales = "free_y") +
# Alternative with facet.grid() (i.e. parameters in rows and objectives in columns) with free y-scales
#facet_grid(as.formula(paste0(facet_var, "~ Objective")), scales = "free_y") +
# Remove x- and y-axis labels
labs(x = NULL, y = NULL, title = plot_title) +
# Choose a light theme
theme_classic() +
# The qualitative 'Set2' palette of RColorBrewer is colorblind friendly
scale_color_brewer(palette = 'Set1')
# Return ggplot object
return(p_out)
}6.3 Plot with Specimen as grouping variable
6.3.1 ISO 25178 height parameters
# 100x objective
# Subset and pivot to longer format for facet plots
Sp_ISO_height_long_100x <- select(sharks_nmp20_100x, all_of(c(x_sp, y_ISO_height, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_height), names_to = "Parameter", values_to = "Value")
str(Sp_ISO_height_long_100x)tibble [525 × 5] (S3: tbl_df/tbl/data.frame)
$ Specimen : chr [1:525] "CC" "CC" "CC" "CC" ...
$ Tooth : chr [1:525] "A" "A" "A" "A" ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 1 1 1 1 1 1 1 1 1 1 ...
$ Parameter: chr [1:525] "Sq" "Ssk" "Sku" "Sp" ...
$ Value : num [1:525] 1.217 -0.332 2.531 2.739 3.419 ...# A tibble: 6 × 5
Specimen Tooth NMP_cat Parameter Value
<chr> <chr> <ord> <chr> <dbl>
1 CC A <10% Sq 1.22
2 CC A <10% Ssk -0.332
3 CC A <10% Sku 2.53
4 CC A <10% Sp 2.74
5 CC A <10% Sv 3.42
6 CC A <10% Sz 6.16 # Define plot
Sp_p_ISO_height_100x <- custom_boxplot(dat = Sp_ISO_height_long_100x, x_axis = x_sp,
plot_title = "ISO height parameters - 100x")
# Print plot
print(Sp_p_ISO_height_100x)Warning: Using shapes for an ordinal variable is not advised
#Same with 20x objective
Sp_p_ISO_height_20x <- select(sharks_nmp20_20x, all_of(c(x_sp, y_ISO_height, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_height), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "ISO height parameters - 20x")
print(Sp_p_ISO_height_20x)Warning: Using shapes for an ordinal variable is not advised
6.3.2 ISO 25178 volume parameters
Sp_p_ISO_vol_100x <- select(sharks_nmp20_100x, all_of(c(x_sp, y_ISO_vol, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_vol), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "ISO volume parameters - 100x")
print(Sp_p_ISO_vol_100x)Warning: Using shapes for an ordinal variable is not advised
Sp_p_ISO_vol_20x <- select(sharks_nmp20_20x, all_of(c(x_sp, y_ISO_vol, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_vol), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "ISO volume parameters - 20x")
print(Sp_p_ISO_vol_20x)Warning: Using shapes for an ordinal variable is not advised
6.3.3 Other ISO 25178 parameters
Sp_p_ISO_others_100x <- select(sharks_nmp20_100x, all_of(c(x_sp, y_ISO_others, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_others), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "Other ISO parameters - 100x")
print(Sp_p_ISO_others_100x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 25 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 25 rows containing missing values (`geom_point()`).
Sp_p_ISO_others_20x <- select(sharks_nmp20_20x, all_of(c(x_sp, y_ISO_others, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_others), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "Other ISO parameters - 20x")
print(Sp_p_ISO_others_20x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 20 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 20 rows containing missing values (`geom_point()`).
6.3.4 Furrow, direction and isotropy parameters
Sp_p_furrow_diriso_100x <- select(sharks_nmp20_100x, all_of(c(x_sp, y_furrow_diriso, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_furrow_diriso), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp,
plot_title = "Furrow, direction and isotropy parameters - 100x")
print(Sp_p_furrow_diriso_100x)Warning: Using shapes for an ordinal variable is not advised
Sp_p_furrow_diriso_20x <- select(sharks_nmp20_20x, all_of(c(x_sp, y_furrow_diriso, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_furrow_diriso), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp,
plot_title = "Furrow, direction and isotropy parameters - 20x")
print(Sp_p_furrow_diriso_20x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 9 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 9 rows containing missing values (`geom_point()`).
6.3.5 SSFA parameters
Sp_p_ssfa_100x <- select(sharks_nmp20_100x, all_of(c(x_sp, y_SSFA, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_SSFA), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "SSFA parameters - 100x")
print(Sp_p_ssfa_100x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 3 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 3 rows containing missing values (`geom_point()`).
Sp_p_ssfa_20x <- select(sharks_nmp20_20x, all_of(c(x_sp, y_SSFA, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_SSFA), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_sp, plot_title = "SSFA parameters - 20x")
print(Sp_p_ssfa_20x)Warning: Using shapes for an ordinal variable is not advised
6.3.6 Save plots
suppressWarnings(
ggexport(plotlist = list(Sp_p_ISO_height_100x, Sp_p_ISO_vol_100x, Sp_p_ISO_others_100x, Sp_p_furrow_diriso_100x,
Sp_p_ssfa_100x),
filename = paste0(dir_out, "/DMTA-Ctenacanths_boxplots-", x_sp, "-100x.pdf"))
)
suppressWarnings(
ggexport(plotlist = list(Sp_p_ISO_height_20x, Sp_p_ISO_vol_20x, Sp_p_ISO_others_20x, Sp_p_furrow_diriso_20x,
Sp_p_ssfa_20x),
filename = paste0(dir_out, "/DMTA-Ctenacanths_boxplots-", x_sp, "-20x.pdf"))
)6.4 Plot with Position as grouping variable
6.4.1 ISO 25178 height parameters
# 100x objective
# Subset and pivot to longer format for facet plots
Pos_ISO_height_long_100x <- select(sharks_nmp20_100x, all_of(c(x_pos, y_ISO_height, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_height), names_to = "Parameter", values_to = "Value")
str(Pos_ISO_height_long_100x)tibble [525 × 5] (S3: tbl_df/tbl/data.frame)
$ Position : chr [1:525] "top" "top" "top" "top" ...
$ Tooth : chr [1:525] "A" "A" "A" "A" ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 1 1 1 1 1 1 1 1 1 1 ...
$ Parameter: chr [1:525] "Sq" "Ssk" "Sku" "Sp" ...
$ Value : num [1:525] 1.217 -0.332 2.531 2.739 3.419 ...# A tibble: 6 × 5
Position Tooth NMP_cat Parameter Value
<chr> <chr> <ord> <chr> <dbl>
1 top A <10% Sq 1.22
2 top A <10% Ssk -0.332
3 top A <10% Sku 2.53
4 top A <10% Sp 2.74
5 top A <10% Sv 3.42
6 top A <10% Sz 6.16 # Define plot
Pos_p_ISO_height_100x <- custom_boxplot(dat = Pos_ISO_height_long_100x, x_axis = x_pos,
plot_title = "ISO height parameters - 100x")
# Print plot
print(Pos_p_ISO_height_100x)Warning: Using shapes for an ordinal variable is not advised
# Same with 20x objective
Pos_p_ISO_height_20x <- select(sharks_nmp20_20x, all_of(c(x_pos, y_ISO_height, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_height), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "ISO height parameters - 20x")
print(Pos_p_ISO_height_20x)Warning: Using shapes for an ordinal variable is not advised
6.4.2 ISO 25178 volume parameters
Pos_p_ISO_vol_100x <- select(sharks_nmp20_100x, all_of(c(x_pos, y_ISO_vol, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_vol), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "ISO volume parameters - 100x")
print(Pos_p_ISO_vol_100x)Warning: Using shapes for an ordinal variable is not advised
Pos_p_ISO_vol_20x <- select(sharks_nmp20_20x, all_of(c(x_pos, y_ISO_vol, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_vol), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "ISO volume parameters - 20x")
print(Pos_p_ISO_vol_20x)Warning: Using shapes for an ordinal variable is not advised
# Alternative with facet.grid()
#Pos_p_ISO_vol <- select(sharks_nmp20, all_of(c(x_pos, y_ISO_vol, grp_colors, grp_shapes, "Objective"))) %>%
# pivot_longer(all_of(y_ISO_vol), names_to = "Parameter", values_to = "Value") %>%
# custom_boxplot(dat = ., x_axis = x_pos, plot_title = "ISO volume parameters")
#print(Pos_p_ISO_vol)6.4.3 Other ISO 25178 parameters
Pos_p_ISO_others_100x <- select(sharks_nmp20_100x, all_of(c(x_pos, y_ISO_others, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_others), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "Other ISO parameters - 100x")
print(Pos_p_ISO_others_100x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 25 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 25 rows containing missing values (`geom_point()`).
Pos_p_ISO_others_20x <- select(sharks_nmp20_20x, all_of(c(x_pos, y_ISO_others, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_ISO_others), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "Other ISO parameters - 20x")
print(Pos_p_ISO_others_20x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 20 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 20 rows containing missing values (`geom_point()`).
6.4.4 Furrow, direction and isotropy parameters
Pos_p_furrow_diriso_100x <- select(sharks_nmp20_100x, all_of(c(x_pos, y_furrow_diriso, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_furrow_diriso), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos,
plot_title = "Furrow, direction and isotropy parameters - 100x")
print(Pos_p_furrow_diriso_100x)Warning: Using shapes for an ordinal variable is not advised
Pos_p_furrow_diriso_20x <- select(sharks_nmp20_20x, all_of(c(x_pos, y_furrow_diriso, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_furrow_diriso), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos,
plot_title = "Furrow, direction and isotropy parameters - 20x")
print(Pos_p_furrow_diriso_20x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 9 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 9 rows containing missing values (`geom_point()`).
6.4.5 SSFA parameters
Pos_p_ssfa_100x <- select(sharks_nmp20_100x, all_of(c(x_pos, y_SSFA, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_SSFA), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "SSFA parameters - 100x")
print(Pos_p_ssfa_100x)Warning: Using shapes for an ordinal variable is not advisedWarning: Removed 3 rows containing non-finite values (`stat_boxplot()`).Warning: Removed 3 rows containing missing values (`geom_point()`).
Pos_p_ssfa_20x <- select(sharks_nmp20_20x, all_of(c(x_pos, y_SSFA, grp_colors, grp_shapes))) %>%
pivot_longer(all_of(y_SSFA), names_to = "Parameter", values_to = "Value") %>%
custom_boxplot(dat = ., x_axis = x_pos, plot_title = "SSFA parameters - 20x")
print(Pos_p_ssfa_20x)Warning: Using shapes for an ordinal variable is not advised
6.4.6 Save plots
suppressWarnings(
ggexport(plotlist = list(Pos_p_ISO_height_100x, Pos_p_ISO_vol_100x, Pos_p_ISO_others_100x, Pos_p_furrow_diriso_100x,
Pos_p_ssfa_100x),
filename = paste0(dir_out, "/DMTA-Ctenacanths_boxplots-", x_pos, "-100x.pdf"))
)
suppressWarnings(
ggexport(plotlist = list(Pos_p_ISO_height_20x, Pos_p_ISO_vol_20x, Pos_p_ISO_others_20x, Pos_p_furrow_diriso_20x,
Pos_p_ssfa_20x),
filename = paste0(dir_out, "/DMTA-Ctenacanths_boxplots-", x_pos, "-20x.pdf"))
)7 Plot anisotropy vs. complexity
7.1 Facetting around Specimen
# set up plot
Sp_p_bi <- ggplot(sharks_nmp20, aes(x = Asfc, y = epLsar)) +
# Scatterplot
geom_point(mapping = aes(color = .data[[grp_colors]], shape = .data[[grp_shapes]]), size = 4, alpha = 0.7) +
# Adjust axes labels
labs(x = "Complexity (Asfc)", y = "Anisotropy (epLsar)") +
# The qualitative 'Set2' palette of RColorBrewer is colorblind friendly
scale_color_brewer(palette = 'Set1') +
# Grid with Objectives in rows and Specimen in columns
facet_grid(Objective ~ Specimen) +
# Choose a light theme
theme_bw()
# Print plot
print(Sp_p_bi)Warning: Using shapes for an ordinal variable is not advised
7.2 Facetting around Position
# set up plot
Pos_p_bi <- ggplot(sharks_nmp20, aes(x = Asfc, y = epLsar)) +
geom_point(mapping = aes(color = .data[[grp_colors]], shape = .data[[grp_shapes]]), size = 4, alpha = 0.7) +
labs(x = "Complexity (Asfc)", y = "Anisotropy (epLsar)") +
scale_color_brewer(palette = 'Set1') +
facet_grid(Objective ~ Position) +
theme_bw()
# Print plot
print(Pos_p_bi)Warning: Using shapes for an ordinal variable is not advised
8 PCA
8.1 Prepare data and select parameters
# Remove rows with NA (complete cases)
data_pca <- na.omit(sharks_nmp20)
# Convert grouping variables into factor()
data_pca[["Specimen"]] <- factor(data_pca[["Specimen"]])
data_pca[["Tooth"]] <- factor(data_pca[["Tooth"]])
data_pca[["Position"]] <- factor(data_pca[["Position"]])
str(data_pca)'data.frame': 93 obs. of 43 variables:
$ Specimen : Factor w/ 2 levels "CC","MM": 1 1 1 1 1 1 1 1 1 1 ...
$ Tooth : Factor w/ 9 levels "A","B","C","D",..: 1 1 1 1 1 1 1 1 1 2 ...
$ Location : chr "loc1" "loc1" "loc1" "loc2" ...
$ Objective : chr "100x" "100x" "100x" "100x" ...
$ Measurement : chr "meas1" "meas2" "meas3" "meas1" ...
$ Position : Factor w/ 2 levels "bottom","top": 2 2 2 1 1 1 1 1 1 2 ...
$ NMP : num 3.03 3.05 3.34 9.84 9.78 ...
$ NMP_cat : Ord.factor w/ 2 levels "<10%"<"10-20%": 1 1 1 1 1 1 2 2 2 1 ...
$ Sq : num 1.22 1.18 1.16 1.66 1.66 ...
$ Ssk : num -0.332 -0.345 -0.355 1.53 1.454 ...
$ Sku : num 2.53 2.45 2.42 8.22 7.72 ...
$ Sp : num 2.74 2.73 2.69 8.77 8.39 ...
$ Sv : num 3.42 3.25 3.17 3.39 3.45 ...
$ Sz : num 6.16 5.98 5.86 12.16 11.84 ...
$ Sa : num 0.995 0.969 0.956 1.14 1.154 ...
$ Smr : num 5.865 4.991 4.781 0.361 0.405 ...
$ Smc : num 1.5 1.47 1.43 1.3 1.32 ...
$ Sxp : num 2.62 2.5 2.48 2.77 2.76 ...
$ Sal : num 20.7 20.4 20.3 15.3 15.2 ...
$ Str : num 0.557 0.53 0.537 0.408 0.396 ...
$ Std : num 58.5 58.5 58.5 84.5 84.5 ...
$ Ssw : num 0.658 0.658 0.658 0.658 0.658 ...
$ Sdq : num 0.179 0.176 0.176 0.409 0.395 ...
$ Sdr : num 1.51 1.47 1.46 4.86 4.71 ...
$ Vm : num 0.0382 0.033 0.032 0.2068 0.2047 ...
$ Vv : num 1.54 1.5 1.47 1.51 1.53 ...
$ Vmp : num 0.0382 0.033 0.032 0.2068 0.2047 ...
$ Vmc : num 1.21 1.18 1.15 1.16 1.18 ...
$ Vvc : num 1.39 1.35 1.32 1.35 1.37 ...
$ Vvv : num 0.148 0.148 0.147 0.154 0.155 ...
$ Maximum.depth.of.furrows: num 2.33 2.29 2.29 4.51 4.52 ...
$ Mean.depth.of.furrows : num 0.834 0.808 0.792 1.292 1.286 ...
$ Mean.density.of.furrows : num 2027 2022 2034 2211 2229 ...
$ First.direction : num 37 37 37 90 90 ...
$ Second.direction : num 56.5 90 90 84.3 84.3 ...
$ Third.direction : num 90 56.6 56.6 78.7 78.7 ...
$ Texture.isotropy : num 47.3 51.9 58.7 72.8 68 ...
$ epLsar : num 0.00119 0.00113 0.00124 0.00152 0.0018 ...
$ NewEplsar : num 0.0172 0.0172 0.0173 0.0176 0.0177 ...
$ Asfc : num 1.85 1.79 1.82 9.4 7.54 ...
$ Smfc : num 86.2 48.8 52 37.9 62.8 ...
$ HAsfc9 : num 0.249 0.233 0.172 3.11 1.59 ...
$ HAsfc81 : num 0.591 0.582 0.516 4.706 5.056 ...
- attr(*, "comment")= Named chr [1:36] "%" "µm" "<no unit>" "<no unit>" ...
..- attr(*, "names")= chr [1:36] "NMP" "Sq" "Ssk" "Sku" ...
- attr(*, "na.action")= 'omit' Named int [1:51] 4 5 6 13 14 15 16 19 20 21 ...
..- attr(*, "names")= chr [1:51] "4" "5" "6" "13" ...# Separate data from 20x and 100x objectives
data_pca_20x <- filter(data_pca, Objective == "20x")
head(data_pca_20x) Specimen Tooth Location Objective Measurement Position NMP NMP_cat
1 CC B loc2 20x meas1 top 1.548641 <10%
2 CC B loc2 20x meas2 top 1.543759 <10%
3 CC B loc2 20x meas3 top 1.529900 <10%
4 CC C loc3 20x meas1 top 4.550709 <10%
5 CC C loc3 20x meas2 top 4.539939 <10%
6 CC C loc3 20x meas3 top 4.706313 <10%
Sq Ssk Sku Sp Sv Sz Sa Smr
1 3.812817 -0.25826476 3.040479 12.67517 12.24684 24.92201 3.010090 0.05951543
2 3.828599 -0.25000042 3.051720 13.19451 12.25975 25.45426 3.024180 0.06553555
3 3.835455 -0.23749999 3.056678 13.86100 12.12212 25.98312 3.029954 0.03624208
4 3.380547 0.05532171 2.746938 11.10437 13.34968 24.45405 2.757626 0.09431685
5 3.407036 0.07149632 2.688279 10.42104 13.45370 23.87474 2.788465 0.21640865
6 3.403758 0.05670066 2.670761 10.33245 13.52249 23.85494 2.786099 0.15631119
Smc Sxp Sal Str Std Ssw Sdq Sdr
1 4.849032 8.469783 74.41928 0.5539008 100.01447 4.39298 0.2275007 2.395721
2 4.844832 8.504251 76.06079 0.5588543 100.01591 4.39298 0.2275266 2.396175
3 4.832575 8.499293 77.44267 0.5635176 100.01410 4.39298 0.2271957 2.388276
4 4.464514 6.246341 65.32191 0.3280670 94.00976 4.39298 0.2269893 2.432957
5 4.559454 6.196548 63.94042 0.3138926 94.00786 4.39298 0.2310735 2.521945
6 4.568143 6.259208 64.40647 0.3120792 94.00386 4.39298 0.2339082 2.581086
Vm Vv Vmp Vmc Vvc Vvv
1 0.1471205 4.996153 0.1471205 3.311605 4.475210 0.5209426
2 0.1497871 4.994619 0.1497871 3.333089 4.473377 0.5212422
3 0.1529777 4.985553 0.1529777 3.340962 4.465460 0.5200924
4 0.1467542 4.611273 0.1467542 3.118653 4.261227 0.3500466
5 0.1436818 4.703140 0.1436818 3.184156 4.368114 0.3350264
6 0.1410312 4.709175 0.1410312 3.179453 4.372397 0.3367782
Maximum.depth.of.furrows Mean.depth.of.furrows Mean.density.of.furrows
1 14.21726 3.632181 986.8928
2 14.09328 3.631081 985.8758
3 14.30434 3.607165 988.6613
4 16.82153 4.914502 1157.0846
5 17.00578 4.866629 1162.9112
6 17.10022 5.112306 1168.0118
First.direction Second.direction Third.direction Texture.isotropy epLsar
1 114.01782 101.3849 106.9970 18.08598 0.001686338
2 114.01934 101.3801 106.9979 18.18273 0.001683542
3 114.00965 101.3838 106.9971 18.20382 0.001626999
4 90.00014 143.0005 106.9825 27.24073 0.003479588
5 89.99183 142.9990 106.9935 26.31377 0.003550548
6 89.98081 142.9995 101.4261 25.47079 0.003555539
NewEplsar Asfc Smfc HAsfc9 HAsfc81
1 0.01709972 4.213370 9.813006 0.2784421 0.5725004
2 0.01711588 4.181272 9.171428 0.2837376 0.5886052
3 0.01712405 4.158485 9.171428 0.2987777 0.5815200
4 0.01631623 4.912851 12.019806 0.2579861 0.4290443
5 0.01626168 5.137305 10.499466 0.2725659 0.4088598
6 0.01625295 5.296622 9.813006 0.2831087 0.4116543 Specimen Tooth Location Objective Measurement Position NMP NMP_cat
1 CC A loc1 100x meas1 top 3.029439 <10%
2 CC A loc1 100x meas2 top 3.049307 <10%
3 CC A loc1 100x meas3 top 3.344563 <10%
4 CC A loc2 100x meas1 bottom 9.838464 <10%
5 CC A loc2 100x meas2 bottom 9.778688 <10%
6 CC A loc2 100x meas3 bottom 9.676411 <10%
Sq Ssk Sku Sp Sv Sz Sa Smr
1 1.216901 -0.3324891 2.531089 2.739170 3.419274 6.158444 0.9953219 5.8654996
2 1.179617 -0.3453720 2.445014 2.728211 3.249398 5.977609 0.9692665 4.9913841
3 1.159639 -0.3552218 2.419064 2.691938 3.171959 5.863897 0.9563751 4.7814041
4 1.656338 1.5304097 8.222025 8.770563 3.392293 12.162856 1.1397393 0.3605701
5 1.660077 1.4536571 7.720735 8.385795 3.451509 11.837304 1.1538744 0.4047778
6 1.672438 1.3428177 7.123336 8.221306 3.535302 11.756608 1.1745203 0.3790766
Smc Sxp Sal Str Std Ssw Sdq Sdr
1 1.501282 2.616679 20.72644 0.5567234 58.51426 0.6583791 0.1788874 1.506669
2 1.466514 2.496966 20.41661 0.5302128 58.51580 0.6583791 0.1762588 1.467242
3 1.433144 2.478542 20.34239 0.5374210 58.51426 0.6583791 0.1757932 1.461778
4 1.299871 2.773392 15.28407 0.4082060 84.49350 0.6583791 0.4089966 4.864301
5 1.323791 2.764030 15.15797 0.3964592 84.49088 0.6583791 0.3952110 4.711679
6 1.340443 2.785413 15.04554 0.3796597 84.49234 0.6583791 0.4100509 4.876693
Vm Vv Vmp Vmc Vvc Vvv
1 0.03822744 1.539509 0.03822744 1.208085 1.391427 0.1480817
2 0.03300682 1.499521 0.03300682 1.177222 1.351987 0.1475341
3 0.03196542 1.465109 0.03196542 1.152959 1.318489 0.1466205
4 0.20681889 1.506690 0.20681889 1.155887 1.352251 0.1544399
5 0.20466698 1.528458 0.20466698 1.179001 1.373821 0.1546366
6 0.20529636 1.545739 0.20529636 1.218516 1.389702 0.1560377
Maximum.depth.of.furrows Mean.depth.of.furrows Mean.density.of.furrows
1 2.330753 0.8337031 2026.884
2 2.286745 0.8079501 2022.394
3 2.291874 0.7918662 2033.563
4 4.505557 1.2923442 2210.950
5 4.519094 1.2861711 2228.890
6 4.589052 1.3129296 2198.418
First.direction Second.direction Third.direction Texture.isotropy epLsar
1 37.00334 56.54720 89.98706 47.34345 0.001190735
2 37.00374 90.00751 56.55262 51.88475 0.001128928
3 37.00464 89.99660 56.55126 58.70932 0.001244379
4 90.00308 84.27736 78.66287 72.77926 0.001522344
5 89.99409 84.26899 78.65584 68.00086 0.001796346
6 90.00073 84.26877 142.98863 65.80603 0.001234300
NewEplsar Asfc Smfc HAsfc9 HAsfc81
1 0.01718885 1.851713 86.15295 0.2487842 0.5906368
2 0.01723834 1.787470 48.78199 0.2329559 0.5822596
3 0.01734807 1.822078 51.96430 0.1722548 0.5164332
4 0.01759510 9.401217 37.88585 3.1102284 4.7057083
5 0.01773082 7.539994 62.81191 1.5898377 5.0555013
6 0.01778692 7.913469 51.96430 1.2035583 3.91063758.2 Run PCAs
8.2.1 PCA on all surfaces
8.2.2 PCA on Specimen CC
# Filter data with 100x on Specimen CC and run PCA
#pca_cc_100x <- prcomp(data_pca_100x[data_pca_100x$Specimen == "CC", pca_params], scale. = TRUE)
# Same with 20x
#pca_cc_20x <- prcomp(data_pca_20x[data_pca_20x$Specimen == "CC", pca_params], scale. = TRUE)
# Filter data on Specimen CC and 100x objective
data_pca_100x_tooth <- filter(data_pca_100x, Specimen == "CC")
# Run PCA
pca_cc_100x <- prcomp(data_pca_100x_tooth[ , pca_params], scale. = TRUE)
# Same on 20x
data_pca_20x_tooth <- filter(data_pca_20x, Specimen == "CC")
pca_cc_20x <- prcomp(data_pca_20x_tooth[ , pca_params], scale. = TRUE)8.3 Plots
8.3.1 Eigenvalues
8.3.1.1 On all surfaces
# 100x
pca_all_100x_eig <- fviz_eig(pca_all_100x, addlabels = TRUE, ggtheme = theme_classic(),
title = "PCA all surfaces - Eigenvalues (100x)")
print(pca_all_100x_eig)
# 20x
pca_all_20x_eig <- fviz_eig(pca_all_20x, addlabels = TRUE, ggtheme = theme_classic(),
title = "PCA all surfaces - Eigenvalues (20x)")
print(pca_all_20x_eig)
8.3.1.2 On Specimen CC
pca_cc_100x_eig <- fviz_eig(pca_cc_100x, addlabels = TRUE, ggtheme = theme_classic(),
title = "PCA CC - Tooth - Eigenvalues (100x)")
print(pca_cc_100x_eig)
pca_cc_20x_eig <- fviz_eig(pca_cc_20x, addlabels = TRUE, ggtheme = theme_classic(),
title = "PCA CC - Tooth - Eigenvalues (20x)")
print(pca_cc_20x_eig)
8.3.2 Biplots
8.3.2.1 Plotting function
custom_pca_biplot <- function(dat, datpca, pc = c(1, 2), geom.pt = "point", col.pt, mean.pt = FALSE,
col.pal = brewer.pal(3, "Set2")[1:2], pt.size = 3, pt.shape = 19, pt.fill = "white",
elli = TRUE, elli.type = "convex", repel.lab = TRUE,
col.variable = "black", main.title){
# Define plotting
p_out <- fviz_pca_biplot(dat, axes = pc,
geom.ind = geom.pt, col.ind = datpca[[col.pt]], mean.point = mean.pt,
palette = col.pal, pointsize = pt.size, pointshape = pt.shape, fill.ind = pt.fill,
addEllipses = elli, ellipse.type = elli.type,
repel = repel.lab, col.var = col.variable,
legend.title = col.pt, title = main.title)
# Return plotting object
return(p_out)
}8.3.2.2 On all surfaces with grouping from Specimen (CC vs. MM)
# Biplot of PC1&2 on 100x
pca_sp_100x_12 <- custom_pca_biplot(pca_all_100x, datpca = data_pca_100x, pc = c(1, 2), col.pt = "Specimen",
main.title = "PCA all surfaces - Specimen - PC1&2 (100x)")
print(pca_sp_100x_12)
# Biplot of PC1&3 on 100x
pca_sp_100x_13 <- custom_pca_biplot(pca_all_100x, datpca = data_pca_100x, pc = c(1, 3), col.pt = "Specimen",
main.title = "PCA all surfaces - Specimen - PC1&3 (100x)")
print(pca_sp_100x_13)
# Biplot of PC1&2 on 20x
pca_sp_20x_12 <- custom_pca_biplot(pca_all_20x, datpca = data_pca_20x, pc = c(1, 2), col.pt = "Specimen",
main.title = "PCA all surfaces - Specimen - PC1&2 (20x)")
print(pca_sp_20x_12)
# Biplot of PC1&3 on 20x
pca_sp_20x_13 <- custom_pca_biplot(pca_all_20x, datpca = data_pca_20x, pc = c(1, 3), col.pt = "Specimen",
main.title = "PCA all surfaces - Specimen - PC1&3 (20x)")
print(pca_sp_20x_13)
8.3.2.3 On all surfaces with grouping from Position (bottom vs. top)
# Biplot of PC1&2 on 100x
pca_pos_100x_12 <- custom_pca_biplot(pca_all_100x, datpca = data_pca_100x, pc = c(1, 2), col.pt = "Position",
main.title = "PCA all surfaces - Position - PC1&2 (100x)")
print(pca_pos_100x_12)
# Biplot of PC1&3 on 100x
pca_pos_100x_13 <- custom_pca_biplot(pca_all_100x, datpca = data_pca_100x, pc = c(1, 3), col.pt = "Position",
main.title = "PCA all surfaces - Position - PC1&3 (100x)")
print(pca_pos_100x_13)
# Biplot of PC1&2 on 20x
pca_pos_20x_12 <- custom_pca_biplot(pca_all_20x, datpca = data_pca_20x, pc = c(1, 2), col.pt = "Position",
main.title = "PCA all surfaces - Position - PC1&2 (20x)")
print(pca_pos_20x_12)
# Biplot of PC1&3 on 20x
pca_pos_20x_13 <- custom_pca_biplot(pca_all_20x, datpca = data_pca_20x, pc = c(1, 3), col.pt = "Position",
main.title = "PCA all surfaces - Position - PC1&3 (20x)")
print(pca_pos_20x_13)
8.3.2.4 On Specimen CC with grouping from Tooth
# Biplot of PC1&2 on 100x
pca_tooth_100x_12 <- custom_pca_biplot(pca_cc_100x, datpca = data_pca_100x_tooth, pc = c(1, 2), col.pt = "Tooth",
col.pal = brewer.pal(7, "Set2"), pt.shape = 21,
pt.fill = data_pca_100x_tooth$Tooth, main.title = "PCA CC - Tooth - PC1&2 (100x)")
print(pca_tooth_100x_12)
# Biplot of PC1&3 on 100x
pca_tooth_100x_13 <- custom_pca_biplot(pca_cc_100x, datpca = data_pca_100x_tooth, pc = c(1, 3), col.pt = "Tooth",
col.pal = brewer.pal(7, "Set2"), pt.shape = 21,
pt.fill = data_pca_100x_tooth$Tooth, main.title = "PCA CC - Tooth - PC1&3 (100x)")
print(pca_tooth_100x_13)
# Biplot of PC1&2 on 20x
pca_tooth_20x_12 <- custom_pca_biplot(pca_cc_20x, datpca = data_pca_20x_tooth, pc = c(1, 2), col.pt = "Tooth",
col.pal = brewer.pal(7, "Set2"), pt.shape = 21,
pt.fill = data_pca_20x_tooth$Tooth, main.title = "PCA CC - Tooth - PC1&2 (20x)")
print(pca_tooth_20x_12)
# Biplot of PC1&3 on 20x
pca_tooth_20x_13 <- custom_pca_biplot(pca_cc_20x, datpca = data_pca_20x_tooth, pc = c(1, 3), col.pt = "Tooth",
col.pal = brewer.pal(7, "Set2"), pt.shape = 21,
pt.fill = data_pca_20x_tooth$Tooth, main.title = "PCA CC - Tooth - PC1&3 (20x)")
print(pca_tooth_20x_13)
8.3.3 Save plots
# 100x
suppressWarnings(
ggexport(plotlist = list(pca_all_100x_eig, pca_sp_100x_12, pca_sp_100x_13, pca_pos_100x_12, pca_pos_100x_13,
pca_cc_100x_eig, pca_tooth_100x_12, pca_tooth_100x_13),
filename = paste0(dir_out, "/DMTA-Ctenacanths_PCA-100x.pdf"))
)
# 20x
suppressWarnings(
ggexport(plotlist = list(pca_all_20x_eig, pca_sp_20x_12, pca_sp_20x_13, pca_pos_20x_12, pca_pos_20x_13,
pca_cc_20x_eig, pca_tooth_20x_12, pca_tooth_20x_13),
filename = paste0(dir_out, "/DMTA-Ctenacanths_PCA-20x.pdf"))
)9 sessionInfo()
R version 4.3.2 (2023-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8
[2] LC_CTYPE=English_United States.utf8
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.utf8
time zone: Europe/Berlin
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] lubridate_1.9.3 forcats_1.0.0 stringr_1.5.1 dplyr_1.1.4
[5] purrr_1.0.2 readr_2.1.4 tidyr_1.3.0 tibble_3.2.1
[9] tidyverse_2.0.0 rmarkdown_2.25 RColorBrewer_1.1-3 R.utils_2.12.3
[13] R.oo_1.25.0 R.methodsS3_1.8.2 knitr_1.45 grateful_0.2.4
[17] ggpubr_0.6.0 factoextra_1.0.7 ggplot2_3.4.4
loaded via a namespace (and not attached):
[1] sass_0.4.8 utf8_1.2.4 generics_0.1.3 rstatix_0.7.2
[5] stringi_1.8.3 hms_1.1.3 digest_0.6.33 magrittr_2.0.3
[9] timechange_0.2.0 evaluate_0.23 grid_4.3.2 fastmap_1.1.1
[13] rprojroot_2.0.4 jsonlite_1.8.8 ggrepel_0.9.4 backports_1.4.1
[17] fansi_1.0.6 scales_1.3.0 jquerylib_0.1.4 abind_1.4-5
[21] cli_3.6.2 rlang_1.1.2 crayon_1.5.2 munsell_0.5.0
[25] withr_2.5.2 cachem_1.0.8 yaml_2.3.8 tools_4.3.2
[29] tzdb_0.4.0 ggsignif_0.6.4 colorspace_2.1-0 broom_1.0.5
[33] vctrs_0.6.5 R6_2.5.1 lifecycle_1.0.4 car_3.1-2
[37] pkgconfig_2.0.3 pillar_1.9.0 bslib_0.6.1 gtable_0.3.4
[41] glue_1.6.2 Rcpp_1.0.11 highr_0.10 xfun_0.41
[45] tidyselect_1.2.0 rstudioapi_0.15.0 farver_2.1.1 htmltools_0.5.7
[49] labeling_0.4.3 carData_3.0-5 compiler_4.3.2 10 Cite R packages used
| Package | Version | Citation |
|---|---|---|
| base | 4.3.2 | R Core Team (2023) |
| factoextra | 1.0.7 | Kassambara and Mundt (2020) |
| ggpubr | 0.6.0 | Kassambara (2023) |
| grateful | 0.2.4 | Francisco Rodriguez-Sanchez and Connor P. Jackson (2023) |
| knitr | 1.45 | Xie (2014); Xie (2015); Xie (2023) |
| R.utils | 2.12.3 | Bengtsson (2023) |
| RColorBrewer | 1.1.3 | Neuwirth (2022) |
| rmarkdown | 2.25 | Xie, Allaire, and Grolemund (2018); Xie, Dervieux, and Riederer (2020); Allaire et al. (2023) |
| tidyverse | 2.0.0 | Wickham et al. (2019) |
References
Allaire, JJ, Yihui Xie, Christophe Dervieux, Jonathan McPherson, Javier
Luraschi, Kevin Ushey, Aron Atkins, et al. 2023. rmarkdown: Dynamic Documents for r. https://github.com/rstudio/rmarkdown.
Bengtsson, Henrik. 2023. R.utils:
Various Programming Utilities. https://CRAN.R-project.org/package=R.utils.
Francisco Rodriguez-Sanchez, and Connor P. Jackson. 2023. grateful: Facilitate Citation of r Packages.
https://pakillo.github.io/grateful/.
Kassambara, Alboukadel. 2023. ggpubr:
“ggplot2” Based Publication
Ready Plots. https://CRAN.R-project.org/package=ggpubr.
Kassambara, Alboukadel, and Fabian Mundt. 2020. factoextra: Extract and Visualize the Results of
Multivariate Data Analyses. https://CRAN.R-project.org/package=factoextra.
Neuwirth, Erich. 2022. RColorBrewer: ColorBrewer
Palettes. https://CRAN.R-project.org/package=RColorBrewer.
R Core Team. 2023. R: A Language and Environment for
Statistical Computing. Vienna, Austria: R Foundation for
Statistical Computing. https://www.R-project.org/.
Wickham, Hadley, Mara Averick, Jennifer Bryan, Winston Chang, Lucy
D’Agostino McGowan, Romain François, Garrett Grolemund, et al. 2019.
“Welcome to the tidyverse.”
Journal of Open Source Software 4 (43): 1686. https://doi.org/10.21105/joss.01686.
Xie, Yihui. 2014. “knitr: A
Comprehensive Tool for Reproducible Research in R.”
In Implementing Reproducible Computational Research, edited by
Victoria Stodden, Friedrich Leisch, and Roger D. Peng. Chapman;
Hall/CRC.
———. 2015. Dynamic Documents with R and Knitr. 2nd
ed. Boca Raton, Florida: Chapman; Hall/CRC. https://yihui.org/knitr/.
———. 2023. knitr: A General-Purpose
Package for Dynamic Report Generation in r. https://yihui.org/knitr/.
Xie, Yihui, J. J. Allaire, and Garrett Grolemund. 2018. R Markdown:
The Definitive Guide. Boca Raton, Florida: Chapman; Hall/CRC. https://bookdown.org/yihui/rmarkdown.
Xie, Yihui, Christophe Dervieux, and Emily Riederer. 2020. R
Markdown Cookbook. Boca Raton, Florida: Chapman; Hall/CRC. https://bookdown.org/yihui/rmarkdown-cookbook.